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A new class of turbulence model is described for wall bounded, high Reynolds number 
flows. A specific turbulence model is demonstrated, with results for favorable and adverse 
pressure gradient flowflelds. Separation predictions are as good or better than either 
Spalart Almaras or SST models, do not require specification of wall distance, and have 
similar or reduced computational effort compared with these models. 


Introduction 

One difficulty with current one and two-equation 
turbulence models is the inability to account di- 
rectly for non-equilibrium effects such as those encoun- 
tered in large pressure gradients involving separation 
and shockwaves. Current turbulence models such as 
Spalart ’s one-equation model , 5 the classic k — e and 
Wilcox’s k — to 1 two-equation models have been de- 
signed and tuned to accurately predict equilibrium 
flows such as zero-pressure gradient boundary- layers 
and free shear layers. Application in more complex 
flows can be problematical at best. Although there 
have been many attempts to modify or correct basic 
one- and two-equation models, most of these attempts 
have been only marginally successful in predicting 
complex flows. 

More complex models such as Reynolds stress mod- 
els have been investigated extensively, primarily for 
relatively simple flows but also for complex flows. In 
most cases these models give somewhat better predic- 
tions than the simpler one and two equation models, 
but for complex flows they do not perform much bet- 
ter than the simpler models. One theoretical advan- 
tage of Reynolds stress models is that they directly 
account for non-equilibrium effects in the sense that 
the Reynolds stresses do not respond instantaneously 
to changes to the strain rate but more realistically 
lag them in time and/or space. Unfortunately, The 
Reynolds stress models are usually considerably more 
complicated and numerically stiff than the one- and 
two- equation models, and this has prevented their 
wide application for complex flows. 

In this paper we introduce a new three equation 
model designed to account for non-equilibrium effects 
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without invoking the full formalism of the Reynolds 
» stress models. The basic idea is to take a baseline 
two-equation model and to couple it with a third (lag) 
equation to model the non-equilibrium effects for the 
eddy viscosity. The third equation is designed to 
predict the equilibrium eddy viscosity in equilibrium 
flows. 

We show results obtained with a lag model based on 
the Wilcox k-tu model. Applications to four flows are 
given including am essentially incompressible flat plate 
flows, am essentiadly incompressible adverse pressure 
gradient flow with separation , 2 a transonic bump flow 3 
with a shock wave and separation, a three dimensional 
transonic wing flow . 4 Results using the new model are 
compared with results obtained with Spalart s model 3 
and Menter’s k-u> SST model . 6 Results obtained with 
the new model show encouraging improvements over 
results obtained with the other models. 


The Lag Model 

The differential equations defining the lag model are 


= Pk-eic + Vlli + cTkP'VtVlc) (!) 

= P^-Ccu+Vtp+CTa-P^tVtu) (2) 

p^- = a(R-r) pen ("Vt E ~' v t) ( 3 ) 

where: 


v tF = k/tu 
P k = pv t S 2 
e k = 6*pkw 


R t = Pk/puj 
V w = apS 2 
e w = Ppcu 2 
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a(Rt) - 

^1(Rt + Rt»)J 

(4) 

with parameters 
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The model 

equations 

are composed of an 

underlying 


model, (k — to) with the lag equation augmenting 
the system. The k — to model, given by the first two 
equations is unchanged from the standard model, 
except that is now given by a field PDB, Eq.(3), 
instead of k/cu. 

The boundary conditions axe those of a converted 
scalar: specified on inflow boundaries and extrapo- 
lated on outflow boundaries. At inflow boundaries, k, 
co and v t are set to constants (hereafter referred to as 
koo, Woo, and y too = In this paper, WU 2 * , 

was chosen to be 0.0001, corresponding to a turbulent 
intensity of 0.8%. The value of k external to the wall 
bounded flows is, of course, substantially lower than 
this, since k decays in the absence of mean strains. 
co M is chosen to yield a low value of v tE (here chosen 
as one tenth the molecular v), and *v t is set to v*, . The 
model's predictions are insensitive to the values chosen 
for these constants as is shown in the results. At the 
solid walls a "no lag” boundary condition is enforced: 
the eddy viscosity is set to its equilibrium value: zero 
for hydraulicly smooth walls, finite for rough walls. 
The cases discussed in this paper all use the smooth 
wall conditions. 7 Rough wall boundary conditions are 
also possible in the manner of the k — co model, 2 but 
are not discussed in this paper. 

The turbulent eddy viscosity is governed by Eq.(3). 
This equation simply says that the eddy viscosity goes 
to its equilibrium value (*v tl ) along a streamline like 
a first order dynamical system with a time constant 
given by 1/(au>). The stability of the turbulence 
model as a dynamical system is ensured by the form 
of this source term, given a stable underlying model. 
There is no diffusion term in this equation, and evo- 
lution of the eddy viscosity is dependent only on its 
upstream history and the underlying equilibrium eddy 
viscosity at that point. 

The a term of Eq.3 of the source for the vt Eq.(4) 
governs the amount of lag present in the model. This 
term is made up of three factors. The leading con- 
stant, ao, controls the amount of lag in the model. 
The higher the value of ao, the less lag(shorter time 
constant) in the system, and the closer will follow the 
underlying turbulence model. The second factor in 
Eq.4 goes monotonically from 100 to one as Rt goes 
from 0 to infinity (Fig. 1). This in effect stiffens the 



Fig. 1 Lag Source Rt Factor 

model at low values of Rt so that it will have less lag, 
and will art more like the unmodified model in these 
conditions. In effect, this term "turns on” the lag only 
in turbulent regions (Rt » 1). 

The model requires the storage of one additional 
field variable over the SST model, but does not require 
the wall distance information used in both the SST 
and SA models, thus freeing up the storage required 
for that variable. As the model is computationally 
simple, not requiring the V>t • V>t or Vk • Vto terms 
present in the other models, it actually requires less 
CPU time per iteration than either SST or SA models, 
and similar or fewer iterations for convergence. 

Numerical Method 

The solutions obtained in this paper were done with 
a modified version of the OVERFLOW 8, • code. Fbr 
the mean flow equations, the existing 3rd order upwind 
scheme or the central/matrix dissipation 10 method 
was used. The full Navier Stokes equations were 
solved, as opposed to utilizing the thin layer approxi- 
mation. Converged solutions were indistinguishable in 
terms of skin friction and velocity profiles. Matrix dis- 
sipation was used with 2nd and 4th order smoothing 
coefficients of 2 and 0.1, respectively. The eigenvalue 
limiters were set to zero, and Roe averaging was used 
to form the matrix. Multigrid was employed, both as 
grid sequencing for startup, and during the relaxation 
process. 3 levels of multigrid were used for all solutions 
presented in this paper. 

The turbulence model equations were spatially dis- 
cretized using a 2nd order upwind method with a 
min-mod limiter. This is a departure from the 2 equa- 
tion models implemented in OVERFLOW, which use 
1st order upwind. 1st order upwind was the initial 
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presented in this paper. 

The turbulence model equations was spatially dis- 
cretized using a 2nd order upwind method with a 
min-mod limiter. This is a departure from the 2 equa- 
tion models implemented in OVERFLOW, which use 
1st order upwind. 1st order upwind was the initial 
implementation, but this proved too dissipative, and 
led to excessive grid density requirements for grid in- 
dependent solutions. 

The reason for this can be seen in the 3rd equa- 
tion, which implements the history effects (lag) of the 
model. Using a 1st order upwind on this equation is 
analogous to using a 1st order time integrator to in- 
tegrate an ODE, with grid spacing analogous to the 
ODE time step size. When the 1st order upwinding 
was replaced with a 2nd upwinding, the gridding re- 
quirements dropped back to what would normally be 
required for grid independent solutions with other tur- 
bulence models. The additional work required to solve 
the third equation is offset by the relative simplicity 
of the underlying model. Convergence is rapid and 
robust, as implemented in OVERFLOW. 

Results 

Dissipation of Isotropic Turbulence 

Isotropic turbulence has no mean strain, so that the 
decay of k and tu follow those of the underlying model, 
here the Wilcox k- tu model. The 'Vt equation uncou- 
ples from the other two equations, and the equation 
governing the evolution of the eddy viscosity becomes: 


d^ t 

17 


= a 0 (k - urv t ) 


This decoupling is aesthetically beautiful, from the 
viewpoint of as a ratio of turbulent stess to mean 
strain. When the mean strain vanishes the turbulent 
stresses are not effected by the value of the eddy vis- 
cosity. Similarly, the lag equation does not affect the 
decay of turbulent kinetic energy built into the under- 
lying model under conditions of zero mean strain. 

Equilibrium Channel Flows 

For fully developed channel shear flows, such as Cou- 
ette, or fully developed pipe flows, the model again 
decouples and reproduces the results of the underly- 
ing model. The differential equation becomes: 


Subsonic Flat Plate 

The fine grid for this case is 
101 (streamwise) x 101 (wall normal). Nearly iden- 
tical results were obtained on a 51x51 grid. The wall 
normal grid stretching for these cases was 1.1 and 1.2 
respectively, with initial + spacings of 0.1 and 0.2. 
The initial 4 wall normal points were equispaced. 



Fig. 2 Subsonic Flat Plate Skin Friction 
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Fig. 3 Law of the Wall Velocity Profiles 


u — = a 0 (k-urv t ) 

3x 

and ^ = Q for fully developed Couette or chan- 

nel flows, this simply enforces = k/tu. In the same 
manner as in the decay of isotropic turbulence, if the 
underlying model does a good job under these con- 
ditions (which k - cu does) then the Lag model will 
also. 


For flat plates, the model roughly reproduces the 
original model performance, since the flow is slowly 
varying in the streamwise direction, and the model 
forces *v t to its equilibrium values. Comparison for 
low speed flows (actually M.^ = 0.2, Fig. ) shows a 
good comparison with Karman-Schoener correlation. 
This was expected from the underlying k — cu model’s 
predictions for smooth flat plates. The law of the wall 
is also well reproduced (Fig. 3). 
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wall distance information, and are under consideration 
for use as future baseline models. The y+ require- 
ments were also investigated, and similar behaviour to 
the underlying k — ui model were found. 7 Improve- 
ments and extensions of the which include roughness 
are available for the k— cu model, 1 * and are also under 
consideration for future model improvement. 

Driver CSO Separated Case 

This is a low speed separated case. 3 The axisym- 
metric geometry (Fig. 5) is defined by an external 
streamline determined from experimental data, and 
wall pressures are available in addition to velocity pro- 
files and skin friction. 



Fig. ft Driver CSO Flowfleld 


Upstream boundary conditions were constant total 
pressure and temperature, with static pressure allowed 
to vary and velocity direction aligned with the cylin- 
der axis. The outer streamline was treated as an 
inviscid wall. The viscous wall(the surface of a cylin- 
der) Is a no-slip, adiabatic wall. The downstream 
static pressure was adjusted to match the experimen- 
tal static pressure (Fig. 6) upstream of the interaction 
region(x « -.438m). Upstream length of the cylin- 
drical body was adjusted so that the computed the 
boundary layer thickness at x = -0.438m matched 
the experiment. This “entry length” was the same for 
all three models. 

The standard OVERFLOW low Mach number pre- 
conditioning wcs employed. The grid used in these cal- 
culations was 200(axial) by 160(radial), an e xtremely 
fine grid. The calculations were repeated with a 




100 x 80 grid and demonstrated grid independence in 
the same manner as Bardina et. al. T Surface pressures, 
skin friction and velocity profiles agreed with the fine 
grid results. This case was one which demonstrated 
the need for handling the turbulence model convection 
operator with 2nd order upwind (minmod limiter), as 
the solutions on the two grids are effectively identical 
when the 2nd order operator is employed, but show 
slight differences when computed with a first order up- 
wind operator for the turbulence model. 

The pressure variation predicted by the Lag model 
is closer to the experimental data than either the 
Spalart-Almaras or SST models, although both give 
reasonable agreement. This improvement in pressure 
variation prediction could be important for internal 
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low speed flows with small separations. The model’s 
ability to more accurately predict the pressure varia- 
tions will be repeated in later test eases. 

The skin friction prediction of all four models is 
shown in Fig. 7. All four predict the skin friction rea- 
sonably well, with the SA model predicting a reattach- 
ment point slightly downstream of the experimental 
data and the k — tu model failing to predict separa- 
tion. All models predict Cf too low Just upstream of 
separation, from x = — .1 to x = 0. 

The velocity profiles (Fig. 8) also show good agree- 
ment. The differences between the underlying k — tu 
model and the lag model are relatively small upstream 
of separation, but the lag model predicts the sepa- 
rated profiles more accurately in the separated region 
(x > 0). 

The shear stress profiles are similarly well predicted 
(Fig. 9), for both the evolution of the maximum shear 
stresss and its location. The lag of the model is evident 
most clearly in this figure. Upstream of the separa- 
tion, the shear stress predicted by the lag model is 
true to it’s name, and lags the underlying k— u» model 
especially evidently at the £, = —0.076 station. By 
the i, — 0.101 station, all of the models are predicting 
roughly the same shear stress, even though the velocity 
profile predictions (Fig. 8) show the greatest scatter at 
this location. One of the model’s slight imperfections 
can be seen in Fig.8 and Fig.9, as the outer edge of 
the boundary layer has a kink not seen in the other 
models. 

Bachalo- Johns on Bump 

This test case features the transonic Interaction of 
a fully developed turbulent boundary layer with the 
pressure field created by a circular arc bump on the 
surface of a cylinder. The surface pressure distribu- 
tions for various freestream Mach numbers are avail- 
able, and velocity and Reynolds stress profiles are 
available for the Moo *■ 0.875 case. 

Upstream boundary conditions were constant total 
pressure and temperature, with static pressure allowed 
to vary and velocity direction aligned with the cylin- 
der axis. The outer edge of the flowfield was treated 
by extending the grid 8 bump chords away from the 
wall, and utilizing characteristic(no reflection) bound- 
ary Conditions. The viscous wall(the surface of a 
cylinder) is a no-slip, adiabatic wall. The downstream 
static pressure was held at p^, Upstream length of the 
cylinder was adjusted to match computed boundary 
layer thickness at x = —0.25m to experiment as done 
in the Driver case, and again this “entry length” was 
the same for all three models. 

The grid for this case is 181(streamwise)x78(wall 
normal). Fine grid solutions with a 358 x 161 grid 
were indistinguishable to plotting accuracy, in terms 
of both surface pressures and velocity profiles at both 
Moo — 0.875 and Moo — 0.925 cases. The wall normal 


grid spacing for both normal and fine grids had a y + 
less than 0.17 upstream of the shock. 



Fig. 10 Moo ” 0.875 Skin Friction Insensitivity to 
koo and Rtoo values 


The insensitivity of the solution to freestream 
choices of and koo is illustrated in Fig. 10. Here, 
the predicted x component of skin friction is shown 
for a range of choices of these parameters. The koo 
range corresponds to an initial freestream turbulence 
intensities from- 0.08% to 2.5%. The Rtoo range corre- 
sponds to Initial eddy viscosities of from 0.1 to 10~ 3 
molecular viscosity. The surface pressure prediction 
variations are just as insensitive to these variations in 
freestream turbulence values. 

The Lag model reproduces the experimental pres- 
sure distributions (Fig. 11) as well as either the 
Spalart-AIlmaras or Menter SST models, a distinct 
improvement over the underlying k— cu model, which 
consistently misses the shock location and underpre- 
dicts the extent of the flow separation. 

The velocity profiles show the progression of the 
flowfield through the separation(Fig. 12). The Lag 
model predicts a separation point intermediate be- 
tween the predictions of the SST and SA models, and 
has a flow recovery better than the SST, though it still 
does not recover as rapidly as experiment. In these ve- 
locity profiles, there is no obvious ldnk at the edge of 
the boundary layer in contrast to the CS0 flowfield. 

The shear stress profiles(Fig. 13) show the Lag 
model’s increased prediction of T max downstream of 
the separation point, though it is not as large as mea- 
sured in the experiment. Note that this plot has the 
wall distance logarithmic, expanding the inner region 
of the boundary layer. There is a consistent under- 
prediction of the stresses in the inner layer by all of 
the models in the separated region, and all the models 
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tween the predictions of the SST and SA models, and 
has a flow recovery better than the SST, though it still 
does not recover as rapidly as experiment. In these ve- 
locity profiles, there is no obvious kink at the edge of 
the boundary layer in contrast to the CSO flowfield. 

The shear stress profiles(Fig. 13) show the Lag 
model’s increased prediction of T max downstream of 
the separation point, though it is not as large as mea- 
sured in the experiment. Note that this plot has the 
wall distance logarithmic, expanding the inner region 
of the boundary layer. There is a consistent under- 
prediction of the stresses in the inner layer by all of 
the models in the separated region, and all the models 
seriously underpredict the measured maximum shear 
stress. 

ONERA M6 Wing 

The ONERA M6 4 is a venerable 3D test case. In the 
conditions from a = 3° to a = 5° range from a nearly 
attached flowfield to a relatively extensive separation, 
as shown in Fig. 14. 

The grid used in these computations 

was the one shipped as a test case with 
OVERFLOW, which has grid dimensions of 

269(streamwise) x 35(spanwise) x 67 (wall normal) , 

with 201 points streamwise along the wing surface. 
The of the first point off the surface was below 
1.25 over the entire wing surface. No grid resolution 
study was performed. 

The C p predictions of the model are compared with 
experiment in Fig. 15. The Lag model produces a good 
prediction over this entire range of conditions. The 
largest discrepancies in this range are at the wing root, 
and are more likely due to the tunnel wall interference, 
as the flow is attached in this region. As can be seen in 
the ” oil flow” pictures of Fig. 14, all these cases have 
some separation, from a ’’incipient separation” at 3° 
to the rather extensively separated 5° case. 

All four models give good predictions at a = 3° and 
a = 4° conditions, but at a = 5° the separation there 
is an appreciable difference in the predictions provided 
by the various models. The wing tip separation pro- 
gression in particular appears to be well captured by 
both the Spalart-Almaras and the Lag model. The 
SST model has more extensive separation than exper- 
iment, and the separation predicted by the k— tu model 
is less extensive than experiment. 

Discussion 

The lag equation could be coupled to virtually any 
model. We have chosen to couple it to the k — cu 
model for this example. Another possible implementa- 
tion would be to lag the Reynold’s stresses, as opposed 
to the eddy viscosity, via an equation of the form 

Dt 

p" = o/(Rt ) CU (2pt t s ij — T ij) 
to account for anisotropic effects seen in 3D flows. 


The main feature of this new class of models is to 
introduce a lag into the response of the eddy viscosity 
to rapid changes in the mean flowfields so as to em- 
ulate the responses seen experimentally. Virtually all 
turbulence models generate Reynolds stresses that re- 
spond too rapidly to changes in mean flow conditions. 
Even the Reynolds stress models predict overly rapid 
response of the Reynolds stresses to changes in mean 
flow conditions. This is in large part due to the mod- 
els’ need to accurately reproduce equilibrium flows. 

The lag equation gives the existing models an addi- 
tional degree of freedom, without tampering with their 
typically good ability to predict equilibrium flows. 

Summary 

The Lag model gives good results for mild to mod- 
erate 2D separations, and agrees well with 3D cases 
tested. It works well for skin friction prediction at in- 
compressible and supersonic Mach numbers, and pre- 
dicts separation well for incompressible and transonic 
test cases. The model does not require wall distance, 
in contrast to both SST and SA models, and its sim- 
plicity is such that the computational effort is roughly 
equivalent to the simpler 1 and 2 equation models. 
Future work will include efforts to remove freestream 
effects, and will look at free shear flows, along with 
other experimental test cases. 
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a)««3° b)«-4* c)*-5 # 

Fig. 14 Predicted Surface Flow (Lag Model), M » 0.84, Re = 18 x 10* 



Fig. 15 ONERA M6 C p Comparisons, a * 3,4 and 5° 
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